﻿using MathNet.Numerics;
using MathNet.Numerics.Providers.FourierTransform;
using System;
using System.Numerics;
using System.Runtime.CompilerServices;
using UtilZ.Dotnet.Ex.Model;

namespace UtilZ.Dotnet.Ex.FourierTransform
{
    /// <summary>
    /// Fourier变换类
    /// </summary>
    unsafe public partial class FourierEx
    {

        private static void Inverse(in SpanZ<Complex> spectrum, FourierTransformScaling scaling)
        {
            if (spectrum.Length.IsPowerOfTwo())
            {
                if (spectrum.Length >= 1024)
                {
                    Radix2InverseParallel(spectrum);
                }
                else
                {
                    Radix2Inverse(spectrum);
                }
            }
            else
            {
                BluesteinInverse(spectrum);
            }

            switch (scaling)
            {
                case FourierTransformScaling.SymmetricScaling:
                    HalfRescale(spectrum);
                    break;
                case FourierTransformScaling.BackwardScaling:
                    FullRescale(spectrum);
                    break;
            }
        }




        /// <summary>
        /// Bluestein generic FFT for arbitrary sized sample vectors.
        /// </summary>
        private static void BluesteinInverse(in SpanZ<Complex> spectrum)
        {
            SwapRealImaginary(spectrum);
            BluesteinConvolutionParallel(spectrum);
            SwapRealImaginary(spectrum);
        }


        /// <summary>
        /// Swap the real and imaginary parts of each sample.
        /// </summary>
        /// <param name="samples">Sample Vector.</param>
        private static void SwapRealImaginary(in SpanZ<Complex> samples)
        {
            //Complex* ptr = samples.DataPtr;
            //Parallel.For(0, samples.Length, (i) =>
            //{
            //    ptr[i] = new Complex(ptr[i].Imaginary, ptr[i].Real);
            //});
            for (int i = 0; i < samples.Length; i++)
            {
                samples[i] = new Complex(samples[i].Imaginary, samples[i].Real);
            }
        }



        /// <summary>
        /// Radix-2 generic FFT for power-of-two sized sample vectors.
        /// </summary>
        private static void Radix2Inverse(in SpanZ<Complex> data)
        {
            Radix2Reorder(data);
            for (var levelSize = 1; levelSize < data.Length; levelSize *= 2)
            {
                for (var k = 0; k < levelSize; k++)
                {
                    Radix2Step(data, 1, levelSize, k);
                }
            }
        }


        /// <summary>
        /// Radix-2 Step Helper Method
        /// </summary>
        /// <param name="samples">Sample vector.</param>
        /// <param name="exponentSign">Fourier series exponent sign.</param>
        /// <param name="levelSize">Level Group Size.</param>
        /// <param name="k">Index inside of the level.</param>
#if !NET4_0
        [MethodImpl(MethodImplOptions.AggressiveInlining)]
#endif
        private static void Radix2Step(in SpanZ<Complex> samples, int exponentSign, int levelSize, int k)
        {
            // Twiddle Factor
            var exponent = (exponentSign * k) * MathNet.Numerics.Constants.Pi / levelSize;
            var w = new Complex(Math.Cos(exponent), Math.Sin(exponent));

            var step = levelSize << 1;
            for (var i = k; i < samples.Length; i += step)
            {
                var ai = samples[i];
                var t = w * samples[i + levelSize];
                samples[i] = ai + t;
                samples[i + levelSize] = ai - t;
            }
        }

    }


}
